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Abstract 

We consider a ID lattice gas model in which the atoms interact via an infinite number of cluster 
interactions within contiguous atomic chains plus the next nearest neighbor pairwise interaction. 
All interactions are of arbitrary strength. An analytical expression for the size distribution of atomic 
chain lengths is obtained in the framework of the canonical ensemble formalism. Application of 
the exact solution to the problems of self-assembly and self-organization is briefly discussed. 

PACS numbers: 68.65.-k, 68.66.La, 81.16.Dn 
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I. INTRODUCTION 



In recent years there has been growing interest in studies of monatomic chains obtained 
in the processes of heteroepitaxial growth [1], |2], |3| at various substrates. A major goal of 
these studies is to obtain the so-called quantum wires for potential use in microelectronic 
applications. In these experiments one usually aims at obtaining wires of infinite length 
for subsequent use in the studies of the Luttinger liquid ||. In practice, however, one 
frequently encounters the monatomic wires of finite length. These may be useful in more 
practical applications, such as the microelectronic circuitry M. Such finite ID clusters were 
observed, e. g., in Refs. M, H]. Current theoretical interest is in magnetic properties of finite 
monatomic clusters |7| in view of their potential use in magnetic memory devices. The latter 
application, however, will require the development of the techniques of mass production of 
such objects. In this context it would be interesting to study the possibility of their self- 
assembly and self-organization similar to analogous processes in 2D heteroepitaxial systems 

II- 

In the present paper we consider an analytic solution for the cluster size distribution in 
the framework of a ID lattice gas model with an arbitrary number of cluster interactions 
within contiguous atomic chains plus the next nearest neighbor pair interaction. Such model 



can be justified in the framework of the Frenkel-Kontorova model of strained epitaxy (12 
but also can be useful in other cases, e. g., in ab initio approaches, where the cluster 
interactions appear because of the many-body nature of the electron interactions which 
cannot be reduced to the pair interatomic potential. 

II. THE MODEL 

Let us consider a lattice gas model defined by the function E{1) describing the dependence 
of the energy of ID atomic clusters (or chains) on their length I. Such a model was already 



used for the description of the self-assembly phenomena in Ref. [10]. We note, however, 
that in 2D this approach can be only phenomenological because the cluster size does not 
characterize a 2D structure uniquely. In contrast, in ID the cluster is defined unambiguously 
by its length and the function E(l) is also uniquely defined provided the clusters do not 
interact. 
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To apply standard tools of statistical physics we need to find the Hamiltonian correspond- 
ing to our model. To this end we consider the following trial expression 

oo 

H = Vi n i n i+i ■ ■ ■ (i) 

i;l=2 

where V\ are numerical coefficients to be fitted to reproduce the chain energies E(l) and 
rii the occupation numbers taking values or 1 depending on whether site i is empty or 
is occupied by an atom. Assuming we have already fitted VJ_i, VJ_2, . . . , V 2 to the energies 
E{1 — 1) = £//_2, . . . ,E 2 let us consider the chain of length /. From Eq. (HD one finds 

E{1) = E(l — 1)+Vi + Vi- X + + ■ ■ ■ + V 2 , (2) 

where the V terms on the right hand side (r. h. s.) account for the interactions of the newly 
added atom with the rest. Similarly 

E(l - 1) = E(l - 2) + Vi-! + V^ 2 + ...+ V 2 . (3) 

Subtracting Eq. (|]) from Eq. (fj) we get 

E(l) - E(l - 1) = E(l - 1) + Vi - E(l - 2) (4) 

from where it follows that the Hamiltonian of our model is Eq. ([!]) with the coefficients given 
by the recursion relation 

Vx = E x - 2^_! + £,_ 2 (5) 

initialized by Eq = E\ — 0. For the system to be well defined in the thermodynamic limit 
the chain energy E\ cannot grow quicker then linearly when / — » oo. This means that the 
the cluster interactions VJ tend to zero at large I because according to Eq. (|J) they are equal 
to the discrete second derivative of E\ with respect to I. 

In the above model atoms interact only when they belong to the same contiguous chain, 
so that separate chains are not coupled to each other. To make the model more realistic we 
allow for the interchain coupling by adding to the Hamiltonian (P the next nearest neighbor 
(NNN) pair interaction term 

oo 

H = V NNN n i n i+2 + V l n i n i+1 ■ ■ ■ ( 6 ) 

i i;l=2 
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III. THE CANONICAL ENSEMBLE SOLUTION 



We consider a ID lattice with a fixed number of atoms N which can occupy I > N 
lattice sites. Our major goal in this paper is to find the equilibrium distribution of the 
atomic chain sizes at finite temperature T in the model defined by Hamiltonian (||). But 
first we will consider the above problem using a simpler Hamiltonian (|l|), i. e., by neglecting 
the NNN interaction which will be accounted for later. We will seek the canonical ensemble 
solution of our problem by generalizing the approach of Ref. ||13|| . In this approach one first 
have to compute the energy of the system in terms of cluster variables mi — the number of 
clusters of length / and their energy E(l). Because the / — > oo limit is not completely trivial, 
we calculate the configurational energy by assuming that all Vi with / > L are equal to zero. 
In this case the total energy can be calculated as 

8 = J^Ernn + E'lN* - m> ]T(Z - l)V h (7) 

1=1 1=2 

where E' L = J2iL 2 Vl * s the first (discrete) derivative of Ei, my is the total number of clusters 
of length exceeding L—2 which below we will call the big clusters, and iV> = \ N — Ya=i Imi) 
is the number of atoms contained in the big clusters. The last two terms on the r. h. s. 
of Eq. (0) can be verified by first checking their correctness for m > = 1 and then noting 
that cutting a big cluster into two ones amounts to cutting I — 1 couplings of type V\. To 
complete the calculation of the total free energy F tot = E ta t — TS it remains to compute 
the entropy S. The number of atomic configurations can be expressed through the cluster 
variables as 



= — — -A — y ~ * | (8) 

m x \m 2 \ . . .m>! \ m> - 1 

In this equality the first multiplier on the r. h. s. counts all inequivalent permutations 
between the clusters with all big clusters being considered as equivalent while the second 
factor accounts for the number of ways to divide N> atoms into m> big clusters including 
their permutations. We find this factor by observing that the division of A^> atoms into 
clusters of sizes exceeding L — 2 is equivalent to dividing n = N > — (L — 2)m > atoms into 



N>-(L- 2)m> - 1 



k = m > clusters. The latter quantity was calculated in Ref. |L3| to be equal to the binomial 
coefficient 

71 — 1 



/~m~ 1 



k- 1 



provided n > k > 1. Otherwise it is equal to zero. The meaning of this formula is simple. 
In a contiguous chain of n atoms there is n — 1 places to cut the chain into pieces, k pieces 
can be obtained with k — 1 cuts. Hence, the number of the possible cuts of n atoms into k 
pieces is equal exactly to the above combinatorial coefficient. To complete the calculation 
of the total number of configurations we have to multiply Q at by the corresponding factor 
responsible for the configurations of vacant sites [IB] . But before doing this we repair Eq. ([7J) 
by accounting for the omitted NNN pair interaction. To this end we note that as long as the 
contiguous chains are concerned, Vnnn interaction amounts simply to renormalization of 
V3. Additionally this term introduces the interchain coupling between the chains separated 
by a single vacant site (a vacancy). Therefore, the term VjViViV&i> where k\ the number 
of vacancies, should be added to the r. h. s. of Eq. fl7|) (here and above we do not pay 



attention to the boundary conditions |13[ because in this paper we are interested only in 
the thermodynamic limit). So when counting the configurations of empty sites we should 
separate the configurations with different numbers of vacancies. This is achieved by using 
Eq. (§) with the total number of empty sites I — N instead of N, the big clusters in this 
case should exceed L — 2 = 1. Denoting the total number of clusters as k from Eq. (0) we 
get 



n k> 

J "i in. 



I-L-k-1 

L k — k\ — 1 



(9) 



k\\(k - ki)\ 

Our final expression for S is obtained from Eqs. (§) and @ with the use of the Stirling 
formula as 



S L ~ 2 

2c\iac — c/lnQ — vlnv — 2(c — v) ln(c — v) 



k B I i=i 

+ (1-9- c) ln(l - B - c) - (1 - 9 - 2c + v) ln(l - 9 - 2c + v) (10) 
-2c> lnc> + [0> — (L — 2)c>] ln[0> — (L — 2)c>] - [9> - (L - l)c>] ln[0> — (L — l)c>], 

where 6 1 = iV// is the total coverage, q = m^/J is the concentration of clusters of length 
I, c — the total cluster concentration, v = k\jl the vacancy concentration, 9 > = 

9 — Y.f=i lei = iV> / /. This expression for reduced entropy together with the reduced energy 
obtained from Eq. (|7]) are sufficient to obtain the reduced free energy as a function of the 
unknown quantities {q}, c>, and v. Minimizing it with respect to these variables we arrive 
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in the limit L — > oo at the following set of equations 

(c - vf = v{l - 9 - 2c + v) exp ( . (12) 

V k B T J 

To find the above limit it is necessary to know the asymptotic behavior of q at large I. The 
limit was taken by assuming the exponential behavior q oc exp(— XI) which is consistent with 
the resulting equation (|TT|). We note that Eq. ( [12] ) is of the second order in v, so one can 
exclude this variable from the first equation to obtain a closed equation for {q}. From the 
above set, however, it is easier to see that if the NNN interaction is negative, than v — > c as 
T — > 0, i. e., the vacancy concentration becomes equal to the concentration of clusters which 
means that the system becomes ordered: the clusters self-organize into chains separated by 
monovacancies. This is further confirmed by the fact that the entropy Eq. (|K|) tends to zero 
at zero temperature provided Vnnn < and the clusters are size calibrated. Indeed, in the 
limit L — > oo only the first two lines of Eq. fllCf) survive 

s L - 2 

- — = 2c In c — q In c\ — v In v — 2(c — v) ln(c — v) 

+ (1-9- c) ln(l - 9 - c) - (1 - 9 - 2c + v) ln(l -9-2c + v), (13) 

where s = S/I is the entropy per site. 

In connection with this expression it is worth noting that Eqs. ([□]) and ( |T2"D can formally 
be obtained in the thermodynamic limit from the variation with respect to q of the expression 
for the free energy density 

oo oo 
1=1 1=1 

where /i is the Lagrange multiplier. This derivation, however, may cause doubts in the 
case of clusters of sizes I = O(N) when the terms /q = miO(N)/N in the above equation 
acquire discrete values and so are not suitable for the variational treatment. The derivation 
presented above is more rigorous and besides can be modified to be applicable to finite 
systems [0, [L6[]. 

At T = and c — v the entropy density (|13D takes the form 

s L " 2 

— = clnc- 22 qlnq. (14) 



k B ^ 



Thus, if the clusters at T = are size calibrated, i. e., if for some value of I = Iq q = c while 
cj = for all other I, then from the above equation it follows that s = 0. For Vnnn > 
the entropy is positive even at T = meaning disordered state. The issue of size calibration 
will be considered in more detail in section M. 



IV. THE ISING MODEL 

To check Eqs. (|TTD and (|1^) we will apply them to the known exactly solvable problem — 
the ID Ising model. As is known, it is equivalent to the lattice gas model with pair inter- 
atomic interaction. This can be easily shown (see ]HJ] Ch. II) by substituting rii = (1 — o"j)/2 
(<7i = ±1 — the Ising spin variable) into the LGM Hamiltonian 

ij ij 

= \ E(^i/ 4 )^- + Yi^/^i 1 - 2 °j) ( 15 ) 

Thus, the Ising model with NN interaction can be solved in the canonical ensemble formalism 
with the use of our formulas. In this section we will obtain the appropriate formulas in closed 
analytical form. Also, we will compare our results obtained in the canonical ensemble with 
more familiar results obtained in the grand canonical formalism via the transfer matrix 



method. The latter solution can be found e. g., in Ch. V of Ref. ||14j| , where the ID Ising 



model Hamiltonian was considered which in our notation can be written as (cf. Eq. ([15]) 

Him = ( V N at/4) ^Vj<7 i+ i - -h^i^i + (16) 
i 1 i 

was considered. The reduced free energy calculated in the grand canonical ensemble was 
found to be 

Fgc 



In e~ K/i cosh h + \je- K l 2 sinh 2 h + e K ' 2 , (17) 



Nk B T 

where K = VMN/ksT and h = h/ksT is the external field (divided by ksT) which fixes the 
magnetization M. The latter is connected to the coverage 9 by the relation 

9 = (1 - M)/2 (18) 

and can be calculated as 

M = dF GC /dh. (19) 
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To obtain the canonical ensemble free energy we have to perform the Legendre transform 
(see Eq. (|T6|)) 

(20) 



r " F GC + h dFGC 



Nk B T " dh ' 
where the partial derivative is calculated from Eq. (|19"D as 



M 



e~ K / 4 sinh h(l + e"^ 4 cosh h/^e~ K / 2 sinh 2 h + e K l 2 ) 



(21) 



e-^/ 4 cosh h + yJe- K / 2 sinh 2 h + e K / 2 

Finally, to obtain the LGM free energy F in the canonical ensemble formalism we have to 
add the statistical average of the last term on the r. h. s. of Eq. (|T5|): 

F F c 



+ (K/A)(1-2M). 



(22) 



Nk B T Nk B T 

In our formalism the free energy F for the NN LGM can be calculated as follows. Ac 
cording to Eq. (|12"1), in the case Vnnn = 



w = c 2 /(l-#). 



Substituting this into Eq. (|TT| ) we get in the case of the Ising model 



(23) 



Cl = (1 - 9 - c) exp[/2 + (fi-K)(l- 1)], 



(24) 



where /i = fx/ksT. From Eq. (|24|) we easily find 

oo 

c = ^Q=(l-^-c)e7(l-e^) 

oo 

e = £/Q=c + c 2 e- JC /(l-e A -*'). 



(25) 



i=i 



These two equations are sufficient to express the unknown quantities c and /x through the 
coverage 9 and the interaction parameter K — the independent variables in the canonical 
ensemble. From Eq. (p5|) we have after some algebra 



In 



26 



— x + 



26 



x ) + x(l — x) 



2x 



1-40(1 -0)z 



where x = (1 — e These expressions can be used to obtain analytic expressions for all 
quantities of interest. 
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In the canonical ensemble case the free energy is 



fl9 - (1 - 9) In 1 + 



-K 



(26) 



As we see, it is quite different from Eq. (P2|). We were unable to compare the two solutions 
analytically and used the numerical procedure. By fixing h and K we computed M with the 
use of Eq. (pip and F/NksT from Eq. (p2|). Then # was calculated according to Eq. ( [3~8|) 
and F/NksT was computed from Eq. (^). Both values coincided within 13-15 significant 
digits depending on how singular the values of F/NksT are at the chosen values of h and 
K. 

Thus we have solved the ID NN LGM in the canonical ensemble approach and expressed 
all quantities of interest in terms of the atomic density (or coverage) 9 which is more natural 
and easier to measure quantity then the chemical potential variable of the grand canonical 
formalism. 



V. THE SIZE CALIBRATION 



According to Ref. ||10|| , the system exhibits the size calibration of the atomic clusters if 
the energy density per atom E(l)/l has a minimum. The size calibration is easiest to see in 
the case V^nn = when the distribution Eq. ( |TTD takes the form 



[1-9 -c) expfiifd - Ei 



Except for the excluded volume prefactor this expression coincides with the formula proposed 
in Ref. |T(| where a detailed discussion of the size calibration issue is given and which fully 



applies to our case provided the pair attraction Vnn is not too strong. In Ref. [12| we have 
shown that in the framework of the Frenkel-Kontorova model the effective chain energy (in 



fact, the free energy — see Ref. [12| for details) has the form 



E t = V NN {l-l)- (Wi-TS l 



(27) 



chains 



where the last two terms represent the relaxation free energy. The length dependence of 
the relaxation energy W\ is governed by the dimensionless parameter a = k s /k p , where 
the spring constants k s and k p are the second derivatives of the potentials which bind the 
atom to the substrate (k s ) and of the pair interatomic interaction (k p ), respectively. In the 
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FIG. 1: The length dependence of the reduced chain energy with the relaxation energy correspond- 
ing to a = 1CP 5 and Vnn = —0.25 in units of k p f 2 (see the text). The inset shows the location of 
the minimum. 

illustrative calculations below the energy unit was chosen to be equal to the energy of the 
unrelaxed misfit strain k p f 2 , where the interatomic spring constant k p was defined earlier 
and / is the misfit as defined in Ref. [Jl2| . The relaxation entropy Si is practically linear in 
I and so according to Eq. ([5]) essentially contributes only into the NN interaction. 

In our calculations we chosen a = 10~ 5 to be quite small in order to visualize all quali- 
tative details of the self-assembly behavior. Physically this would correspond to very weak 
binding to the substrate. In the above energy units Vnn was chosen to be equal to -0.25 
so that there was a minimum in E{f)/l at I « 100 (see Fig.|l|). The calculation of the size 
distribution for this case is shown in Fig. |^. It qualitatively agrees with the Monte Carlo 
simulations of Ref. |TT| (cf. Fig. 2 of that reference). The most notable feature of the above 
calculations is the change of the position of the maximum of the cluster size distribution 



with lowering temperature — the feature not clearly understood from the theory of Ref. pD 
Although the entropic contribution into Ei is rather small due to the small temperatures 
considered, the shallowness of the minimum in E(l)/l makes it possible that the above shift 
of the maximum is due to the shift in the position of the minimum in E\ = W\ — TSi with 
lowering temperature. To check this possibility we repeated calculations with E\ replaced 
by temperature independent part W\. The results of this calculation is shown in Fig. |3|. The 
only qualitative difference with Fig. ^| is that the bimodal distribution is never observed in 
the case of T-independent case. Presumably, this is due to the ID character of our model, 
because in 2D case where also T-independent interactions were used we found a strongly 
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I 

FIG. 2: Equilibrium distribution of atomic chain lengths at the coverage 9 = 0.1 corresponding to 
(from left to right) l/k B T =20, 22, 30, and 60 (in units of k p f 2 ). The pair interaction V NN = -0.25 
and a = 10~ 5 (for explanation of notation see the text). 
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FIG. 3: The same as in Fig. |2| except that the entropic part of E\ is set to zero and the temperatures 
are (from left to right) l/k B T =12.5, 13, 15, 20, and 60 



bimodal distribution ]TS| 



VI. CONCLUSION 



In this paper we considered an analytical solution of a ID lattice gas model with an 
arbitrary number of cluster interatomic interactions within contiguous atomic chains. We 
illustrated the model by applying it to the description of the equilibrium distribution of 
atomic clusters in strained epitaxy in the thermodynamic limit. However, the model can 
be used also in other cases when the multiatom interactions are of the above type, as well 
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as applied to finite systems, like those discussed in Refs. fL3l |TH|. Our analytical solution 



confirms the observation made in Ref. [DJ in the framework of a phenomenological approach 
that the maximum in the equilibrium distribution of the cluster sizes moves to higher values 
with lowering temperature. This can lead to an interesting kinetic phenomenon consisting 
in the saturated ripening when the system is quenched to lower temperature. 
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